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Abstract. Many biochemical processes can successfully be described by dynamical systems al- 
lowing some form of switching when, depending on their initial conditions, solutions of the dynamical 
system end up in different regions of state space (associated with different biochemical functions). 
Switching is often realized by a bistable system (i.e. a dynamical system allowing two stable steady 
state solutions) and, in the majority of cases, bistability is established numerically. In our point of 
view this approach is too restrictive, as, one the one hand, due to predominant parameter uncertainty 
numerical methods are generally difficult to apply to realistic models originating in Systems Biology. 
And on the other hand switching already arises with the occurrence of a saddle type steady state 
(characterized by a Jacobian where exactly one eigenvalue is positive and the remaining eigenvalues 
have negative real part). Consequently we derive conditions based on linear inequalities that allow 
the analytic computation of states and parameters where the Jacobian derived from a mass action 
network has a defective zero eigenvalue so that — under certain genericity conditions - a saddle-node 
bifurcation occurs. Our conditions are applicable to general mass action networks involving at least 
one conservation relation, however, they are only sufficient (as infeasibility of linear inequalities does 
not exclude defective zero eigenvalues). 

1. Introduction. Many biochemical processes can successfully be described by 
dynamical systems allowing some form of switching, where, depending on, for example, 
initial conditions, solutions of the dynamical system end up in different regions of state 
space (associated with different biochemical functions). Often dynamical systems 
admitting bistability (i.e. the existence of two stable steady states) are used for this 
purpose. There is a long tradition of establishing bistability, both experimentally and 
computationally, in areas ranging from signal transduction (see e.g. [2]) to cell cycle 
(see e.g. [B]). 

From our point of view, however, bistability is too strong a requirement, as al- 
ready a saddle type steady state with just one algebraically simple positive eigenvalue 
and all other eigenvalues having negative real part gives rise to the desired switch- 
ing behaviour (with the global stable manifold of the saddle as a switching surface, 
see pjJJ Remark 3.2]). The approach presented here tries to directly establish such 
points and is hence capable of establishing switching that is not necessarily associ- 
ated to bistability. Therefore we expect this approach to be of particular interest for 
researchers working in Systems Biology and other areas of Quantitative Biology. 

In many applications, bistability of a dynamical system has been established 
numerically using bifurcation analysis or simulations that can become arduous tasks 
even for relatively small systems. Moreover parameter uncertainty is a predominant 
issue in Systems Biology: the dynamical systems consist of a large number of states 
and parameters, while measurement data are often very noisy and data points and 
repetitions are usually few. Hence techniques allowing the direct analytic computation 
of parameter vectors where a given system exhibits switching are desirable. 

In developing these techniques we identified two promising approaches: (i) estab- 
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lishing multiple steady states as a mechanism for possible switching and bistability 
[TU1 [TT1 IT21 [T5] and (ii) establishing points where the dynamical system undergoes 
a saddle-node bifurcation so that the global stable manifold of the saddle is acting 
as a switching surface. The first approach is motivated by the so-called Chemical 
Reaction Network Theory developed by Feinberg and co-workers (see [HI [TjJ [18] and 
[T2] [T9] [20]). The second approach is based on the structure of the Jacobian of a 
mass-action network [9]. This approach was successful for a double-phosphorylation 
mechanism where the nullspace of the Jacobian admits a very special representation 
(cf. 0). 

Here we extend these ideas to mass action networks in general (involving at least 
one conservation relation) by making use of a property that is frequently observed 
in dynamical systems originating in Systems Biology: one often faces dynamical sys- 
tems that involve so-called conservation relations confining trajectories to affine linear 
subspaces of state space. 

Thus, the Jacobian of such a system evaluated at an arbitrary point in state 
space has at least as many zero eigenvalues as there are conservation relations. Con- 
sequently, for a saddle-node to occur at a particular point in state space, the Jacobian 
has to have an additional zero eigenvalue at that point. Generically, mass-action sys- 
tems undergo a bifurcation at that point - one can state conditions guaranteeing a 
saddle-node bifurcation (cf. Section [S] or [9]). One can expect that such sufficient 
conditions for a saddle-node bifurcation can be established for mass action networks 
originating in Systems Biology since there are many parameters which can be cho- 
sen as continuation parameters. Hence, such an additional zero eigenvalue frequently 
entails a saddle-node bifurcation and thus switching in a mass action network. 

The main result of our paper are sufficient conditions guaranteeing such an ad- 
ditional zero eigenvalue that take the form of linear inequality systems and are thus 
easy to check. Moreover, our result is constructive in the sense that the solutions 
to one of the inequality systems determine a state and parameter vector where the 
Jacobian has an additional zero eigenvalue and thus fulfills the necessary degeneracy 
condition for a saddle-node. Infeasibility of all inequality systems does not exclude 
additional zero eigenvalues, hence feasibility of at least one inequality system is a 
sufficient condition for an additional eigenvalue. In case the remaining eigenvalues of 
the linearization have negative real parts such feasibility is generically sufficient for a 
saddle-node and the associated bifurcation into a saddle and a node. We verify this 
splitting in our case studies by computing bifurcation diagrams. 

Finally we'd like to point out that our results are in a certain sense complementary 
to those obtained in [27], [T3] [H [15], [29], [3 [H [5]: all these references present 
sufficient conditions for the global injectivity of a dynamical system defined by a 
biochemical reaction network (not necessarily restricted to mass action systems). In 
particular, these conditions exclude switching. More along the line of our work is the 
approach of Mincheva and coworkers [231 12S] • There the tight connection between the 
characteristic polynomial of the Jacobian and the cycles of certain graphs associated 
to the Jacobian are exploited to derive conditions for certain instabilities (e.g. saddle- 
node or Hopf bifurcations). The major difference to our work is that we are not 
working with the characteristic polynomial but rather exploit the fact (reported in [9]) 
that Jordan blocks of size > 2 imply additional zero eigenvalues (and thus candidates 
for, for example, saddle- node bifurcations). 

In the following we briefly describe the organization of the paper and at the same 
time offer the conclusions that can be drawn. In Section [2] we describe dynamical 
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systems defined by mass action networks, recall some results from [9] and characterize 
positive state vectors where the Jacobian has such an additional and thus defective 
zero eigenvalue (Lemma 12.11 Theorem 12. 4p . Those state vectors arise from elements 
of a semialgebraic set that contains only polynomials of degree two or less - regardless 
of the exponents in the polynomial ODE system defined by a mass action network. In 
Section [3j based on a result from Qualitative Matrix Theory ensuring the existence 
of positive null vectors, we present a sufficient condition allowing the computation of 
elements of that semialgebraic set that takes the form of linear inequality systems. 
The solvability of these inequality systems is then sufficient for the existence of an 
additional zero eigenvalue (Theorem 13.21) . In Section [4j finally, we demonstrate the 
applicability of the results presented here by analyzing as a proof of principle two 
competing mass action networks describing the Gl/S transition in the cell cycle of 
budding yeast. These networks were originally presented in |12j and |19) where their 
investigation was based on subnetwork analysis. Both networks are not accessible by 
the results of [9]. 

For the convenience of the reader we provide some additional information in 
four appendices. In Appendix [X] we recall some remarks concerning saddle-node 
bifurcations in mass-action networks that were made earlier in [9] , in Appendix [Bl 
and Appendix [Q we collect the relevant structural information of the Gl/S transition 
networks discussed in Section |U The final Appendix [D] using basic linear algebra, 
discusses some of the assumptions and results in the present work. 

2. Dynamical systems defined by mass action systems. 

To introduce the notation, we use the network depicted in equation (|2.1j) below. 
This network is analysed in Qj5], where multiple steady states are established. 



A + 2S- 



ki 



AS 2 



B + S- 



k 4 



BS 



AS 2 + BS 



C + 3S 



(2.1) 



A 



B 



kio 



c 



Network (|2.1j) consists of n species (n = 6) and with each species we associate a 
variable Xi representing its concentration and the corresponding unit vector of lR e : 
X\ and ei with A, x 2 and e 2 with B, x% and e^ with S, x^ and e^ with AS 2 , x$ and 
e$ with B S and x§ and eg with C . 

The nodes of the network graph are called complexes and with each complex we 
associate the sum of its constituent species. The above network contains m complexes 
(m = 10): The complex will be denoted by the zero vector £ JR 6 and is used to 
encode that the system is open with respect to A, B and C: A and B can enter and 
leave the system while C can only leave the system. As a complex A is associated 
with ei G lR e , B with e 2 , C with e^, A + 2 S with e% + 2 A S 2 with ei, B + S with 
62 + 63, B S with es, A S 2 + B S with e^ + e§ and C + 3 S with ee + 3 63. 
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The network consists of r reactions (r = 10), e.g. A + 2S — > AS 2 , where the 
complex at the tail of the arrow is called educt complex and the complex at the tip 
of the arrow is called product complex. To each reaction is associated a reaction rate 
Viik, x). For mass action systems i»i(fc, x) is proportional to the product of (powers of) 
concentrations of the species forming the educt complexes: let yi be an educt complex 
vector, then one has Vi(k, x) = ki x Vi (where x p = Y[j %j for n- vectors x and p). For 
the above network one obtains 



v(k,x) = (kixixl, k 2 x A , k 3 x 2 x 3 , k t x 5 , k 5 X4,x 5 , k 6 xi, fc 7 , fc 8 , k 9 x 2l k w x e ) T . 

We collect the exponents yi of the monomials contained in Vi(k, x) in the rate- exponent 
matrix y. For the above network one obtains the (n x r)-matrix 

y = [yi, - ,2/io] 

= [e a + 2e 3 , e 4 , e 2 + e 3 , e 5 , e 4 + e 5 , e x , 0, 0, e 2 , e 6 ] . 

The reactions are encoded in the stoichiometric matrix 5 1 , where each column corre- 
sponds to one reaction and is defined as the difference between product and educt com- 



plex. For example for the reaction A + 2 S — > A S 2 one obtains r\ 
The stoichiometric matrix for the above network is 



-(ei + 2e 3 ) + e 4 . 



S = 



1 


2 
-1 





A reaction network then defines a dynamical system 

x — S v(k, x), 

which in case of the above network translates to 



(2.2) 



x\ = kj — k^xi — k\xix 3 + k 2 x^ 

x 2 = kg — k 9 x 2 - k 3 x 2 x 3 + kiX 5 

±3 = —k 3 x 2 x 3 — 2kiXix^ + 2k 2 X4 + /C4X5 + 3^X4X5 

±4 = kiXix 3 — k 2 X4 — ^53:4X5 

x 5 = k 3 x 2 x 3 - fc 4 x 5 - k 5 X4X 5 

x e = k 5 X4X 5 - k w x 6 . 



In general we consider a mass action network with n species, m complexes and r 
reactions. Any such system defines a dynamical system in the form given in (|2.2j) . 
Note that v(k,x) £ M r is a monomial, vector- valued function of the form 

v(k,x) = diag(fc) <fi(x), 

where diag(fc) is a (r x r) diagonal matrix with the ki on the diagonal and 4>(x) = 
i xVi )i=i r ^ ^ r IS a vec tor of monomials in x. Note that the rate-exponent matrix 
y defined above contains the exponent vectors of the monomials contained in <j){x). 
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In the sequel, we speak of steady states (k,x) of (|2.2p when Sv(k,x) vanishes for 
positive (k, x). 

For many realistic systems in Systems Biology the matrix S G JR nxr does not 
have full row rank s := rank(5) (i.e. s < n). This gives rise to n — s conservation 
relations: let Z be any matrix whose columns form a basis of kei(S T ), the left kernel 
of S. Solutions x(t) to (|2.2|) then satisfy 

Z T x{t) = Z T x(0) =: c, (2.3a) 

that is, these solutions lie in invariant domains cc(0)+im(S l ) that are parallel translates 
of im(5). For the above example (12.11) one obtains 

X3 + 2x4 + %5 = c . (2.3b) 

2.1. The Jacobian associated to a mass action network. At positive (k, x) 
the Jacobian of a mass action network (by this we mean the Jacobian of a dynamical 
system defined by a mass action network) is given by: 

Jac{k, x) = S diag (v(k, x)) y T diag (aT 1 ) , (2.4) 

with stoichiometric matrix S and rate-exponent matrix y. 

Observe that a positive pair (k, x) is a steady state of (|2.2j) if and only if v(k, x) G 
int (ker (S) D -R> ) (where int(-) denotes the relative interior). The pointed polyhe- 
dral cone ker (S) PI K> is generated by a finite set of unique (up to scalar multipli- 
cation) extreme rays |26) . The calculation of these rays is in general computationally 
hard, however, there exists a variety of algorithms and software tools, for example 
[2"2l |3"0] . Let p be the number of extreme rays and let E be a matrix whose columns 
are generators of ker (S) Pi -K> . Then (fc, x) is a positive steady state if and only if 
there exists a v with 

v(k,x) =Ev >0, v G J?> . (2.5a) 

So we ask for all components of E v to be (strictly) positive. We collect all such v in 
the set 

V := \v e R p >a \Ev > oj . (2.5b) 

Since E is a nonnegative matrix, V consists of the positive orthant JR^, , i.e. the 
interior of -K> , and potentially certain faces of -K> (i.e. elements v S V are either 
positive or nonnegative with predefined sign pattern). 

As we are interested in the Jacobian Jac(k, x) evaluated at a positive steady state 
we use (|2.5ap in \2A\ to obtain 

Jac{k,x) = J{v,x) = N(v) diag (a;' 1 ) , (v, x) e V x , (2.5c) 

with the z/-linear 

iV(i/) := 5 diag (E 1/) ^ T G M nxn , (1/, as) G V x J^ . (2.5d) 

We'd like to emphasize that points (l>, x) G V x -2?™ define points (fc, x) G -K>o x -K>o 
via 



k = diag(0(a;^ 1 )) £ v 
5 



(2.6) 



Hence finding points {y, x) where J(v, x) is singular is equivalent to finding points 
(k, x) where the Jacobian Jac(k, x) is singular. Null vectors of Jac(k, x) = J{y, x) of 
the form diag (x) z will be obtained from the identity 



J(v, x) diag(x) z = N(v) z = H(z) v = , v e V , 



with the z-linear 



H{z) := SAmg(y T z)E G M 



nxp 



(2.7a) 



(2.7b) 



Our goal in (I2.7a[) is to use a condition from Qualitative Matrix Theory that entails 
the existence of a positive null vector v for the matrix H(z) (cf. Theorem 13. 2p . 

Example 1 (J(y,x) derived from network (|2.1[0 . The generator matrix of 
ker (5*) n 5?> is given by 



E = 



1 











1 


1 

















1 








1 





1 























1 








1 
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1 
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1 











1 

















1 



(2.8a) 



and hence satisfies 
The matrix J{v, x) is given by 



Ev > <^> v > and thus V = M> . 



(2.8b) 



J(u, x) 








2(1/1+1/5) 


"1 










1^2+1/4+1/5 


x ? 











1/2+1/5 





l'2 





2(1/1+1/5) 


X2 


X3 




15 




1/2+1/5 


41/1+1/2 +5l/5 




1/2+31/5 





xi 


X2 




X4 


xs 







2(i/i+ 3 i/ 5 ) 









x\ 


X3 


X4 


Xs 





1/2+1/5 


1/2+1/5 


_E5_ 


1/2+1/5 





X2 


X3 


X4 


X5 











£i 


i& 


— EZl 


X4 


x 5 


x 6 - 



2.2. Zero eigenvalues of the Jacobian of a mass action system. We as- 
sume s — rank(S) < n so that the Jacobian always has n — s zero eigenvalues. In 
addition we assume 

im(S) = ±m(J(v,x)) . (2.9) 

In other terms, we assume the columns of the matrix Z from (|2.3ap to form a basis 
for ker ( J T {v, x)) so that J(z^, x) does not possess more conservation laws than S. In 
the end, we will have to validate this condition (|2 .Q[) (cf. Appendix lD.il and ID.3|) . 

In looking for bifurcations, we reduce the system to the afhne subspaces x(0) + 
im(S'). To this end let U, W be orthonormal bases of im(S), im^S 1 )- 1 -, respectively 
and introduce 

£ = U T x, r) = W T x and x{£,r[) = U£, + Wrj (2.10a) 
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to obtain the reduced system 



i = U T S v (fc, 77)) =: <? (£, 77, fc) (2.10b) 
77 = 0. (2.10c) 

Then the upper left block of the Jacobian of this mass action network is given by 

D s g(£,ri,k) = U T Jac(k,x(C, V ))U 

and at (u, x) G V x J?™ by 

G(£, 77, v) = U T J{y, x(t 77)) U g iR sxs , (2.11) 

where we recall the relation (|2.6p between fc, ^ and a;. In [5] we presented a method 
that links zero eigenvalues of G(£, 77, to zero eigenvalues of J{y, x(£, 77)). Lemma l2.ll 
below is required for Theorem 12.41 the main result of this section. We state it here 
without proof, for a proof see [5] . 

We start with some notation and, as in |31| . call an eigenvalue A of a matrix 
A G M nxn defective if its algebraic multiplicity m a i g (X) is greater than its geometric 
multiplicity m geo (X), that is, if the multiplicity of A as a root of the characteristic poly- 
nomial is greater than the number of linear independent eigenvectors corresponding 
to A. Hence, Ao = is a defective eigenvalue of A if and only if dim(ker(j4) + im(A)) 
is less n. This can be stated in the following way: 

Fact 1. Ao = is a defective eigenvalue of a matrix A £ M nXn iff there exists 
an x =/= with x G im (A) D ker (A). 

Remark 1. An alternative argument for Fact[l]is based on the Jordan Canonical 
Form of a matrix A (cf., for example, J31f). Assume annxn matrix A with eigenvalue 
Ao = and m a i g (Ao) > m geo (\o) in Jordan Canonical Form. Then the m a i g x m a i g 
block matrix corresponding to Xq is not the zero-matrix, implying the existence of 
nontrivial U\ ^ ui with Au\ = and Au2 — u± and hence U\ G ker(yl) (~| im(^4). 

We recall another fact from Lemma 1 in [9]: 

Lemma 2.1. Let A G M nxn be a matrix of rank s < n and let U be orthonormal 
basis for im(A). Then Ao = is a defective eigenvalue of A if and only if Xq =0 is 
an eigenvalue of Bi := U T AU G 1R SXS . 

Based on Fact [1] and Lemma \2. II one is led to following observation: 

Lemma 2.2. Let Zq be a basis ofim(S) ± . Then the Jacobian G{^,rj,v) of the 
reduced system, evaluated at v G V and x = U £ + W 77 G -R™ (cf. H2.ll]) and i 2.10a)) ), 
has a zero eigenvalue if and only if there exist a nontrivial vector z G M n , a vector 
x G JR"o an d a vec t° r v £ V with 

H{z)v = Q (2.12a) 
Zl diag(ai) z = 0. (2.12b) 

In the sequel, we take for Zq the matrix Z describing the conservation laws (cf. 

Proof. From Lemma |2~T1 follows that G(£, 77, v) has Ao = as an eigenvalue, if and 
only if J(v, x) has Ao = as a defective eigenvalue. From Fact [1] follows that J(y, x) 
has a defective eigenvalue, if and only if there is a nontrivial vector z G ker (S')nim (S). 
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That is, S must satisfy N(v) diag (|) z = and Z T z = (cf. (|23c|) and ([23)1 ). Let 
z = diag (x) z, then (I2.12al) and (|2.12b|) follow. □ 

First we consider condition (|2.12bj) and establish necessary and sufficient condi- 
tions for the existence of solutions (x, z) G -K>o x -K™, where we assume that z is 
given. 

Lemma 2.3. Let M G Rnxn be any matrix and let z G M qi be given. Then there 
exists a positive vector x G 1R>q such that 

M T diag(a;) z = 0, 

if and only if 

3a; G ker (M T ) with sign(w) = sign(z). (2.13) 
In this case x — (x) i=1 is given by 



Xi > 0, arbitrary, if Zi = 0. 



(2.14) 



Proof. Assume M T diag(x) 2 = holds for positive x and some z. Then u> := 
diag(x) z G ker (M T ^j and sign(w) = sign(z) follows from positivity of x. Vice versa, 
let z G M and ui G ker (M T ^j with sign(w) = sign(z) be given. Let x be as in (|2.14[) . 
Then sign(cj) = sign(z) implies positivity of x and one has diag(a;) z = lo G ker (Af T ) . 
□ 

Remark 2. Observe that, given a vector z, the condition i2.13\) takes the form 
of linear inequalities: one has to establish feasibility of the system 

M T ijj = 0, sign (z^ uji > 0, if Zi =/= and loi = 0, if Zi = 0. 



Remark 3 (Connection to [9 ). The condition h2.12a\) requires the symbolic 
computation of ker (H (is)) . This can be of forbidding complexity, especially for large 
networks, even though it is in principle possible. 

So far, the only application of the simple fact in Lemma \2.1\ we are aware of was in 
J3[. There we analysed a mass action network describing the double phosphorylation 
of a protein. For this network we obtained a symbolic representation of ker (N(v)) 
that could be brought into a v -independent form. In general, the previous approach 
requires positive solutions to some well-defined polynomial equations in v and is thus 
limited to certain classes of systems (cf. \ 10$). 

In the sequel, we employ the structure of H(z), given by (|2.7bj) . when discussing 
H(z) v = Q. Let the columns of 5 G M rx( - r -^ be a basis of ker(S) and let S# G M rxr 

be a matrix such that S# So = ; . If we let the columns of S c be a basis 

for im(S' T ) and if we denote the Moore-Penrose inverse (Sj Sq)~ 1 S t by Sf we will 

consider a particular such S# by setting S^ rt = 

Equation (|2.12ap is now equivalent to 



5 q = diag(^ T z) Ev 



(2.15) 



for some vector a G IR r s and, by left multiplication with S^ r , to 

= ^ or *diag (y T z)E. 



a 




P(z) 

Q{z) 



v with 



P(z) 

Q{z) 



(2.16) 



Observe that z and v satisfy H(z) v = (cf. ([2.7aj) &; (|2.7b[) ) if and only if one has 

Q(z)v = 0, i/eV. (2.17) 

The corresponding a will be given by P(z)v. We note that the elements of the matrices 
P{z) and Q(z) are linear forms in z. Appendix ID.2I shows that the condition (I2.17[) 
is independent from the chosen bases for ker(S') and im(S T ). 

Theorem 2.4. The Jacobian G(£, n, v) of the reduced system, evaluated at £ and 
r\ with x — x(^,rj) G IRyQ as in \2. 1 0a\) and v G V, has zero as an eigenvalue with 
algebraic multiplicity > 1 if and only if there exist z E IR n , u) € lR n and fi G -K>o 
with 



Q(z)n = 0, E[i > 0, Z t lu = 0, sign(w) = sign(z) 



(2.18) 



Proof. With the settings uj = diag(a;)z as in (|2.14l) . x = U£ 
and v = fx, Theorem 12.41 follows immediately from the Lemma 
the equivalence of H(z)v = with (|2.17[) . □ 



W-n as in (I2.10a.p 
Lemma 12.31 and 



Remark 4 (Open condition (|2.16|1 ). Observe that Q(z) [i = in 112.18]) can also 
be written in the form Q(/i) z = since N(fi)z — H(z)fi ( cf. \2. 7a\ ) with \2. 5d\) & 
\2. 7b\ ) ) implies the equivalence of 12. 1 6\) and 



a 






z with 









. Q(aO . 





(2.19) 



This reformulation reveals that \2. 1 &jl is an open condition: Given a particular solu- 
tion (z,Q,fi), there will exist a solution (z,uj,fi) for all fi's that are sufficiently close 
to fj,. So, there is some freedom in the choice of \i, cf. Avvendix \D.3[ 



Note that the semialgebraic set given by (|2.18p is always defined by polynomials 
of degree two or less, independent of the exponents in the polynomial ODEs. Any 
element gives rise to a defective eigenvalue of the Jacobian Jac(k,x). For the 
computation of elements of that set we will later on employ the following observation: 
in case the vector fi in (|2.18p can be chosen as a positive null vector of Q(z), the 
condition Efi > is automatically satisfied. Thus we arrive at a sufficient condition 
for a defective eigenvalue of the Jacobian Jac{k,x) by imposing conditions on the 
matrix Q(z) that imply the existence of a positive null vector fi and conditions on the 
vector z ensuring the sign-compatibility of z with ker (Z T ). 

Since the elements of Q(z) derived from a mass action network are always linear 
forms in z, one can determine all sign patterns that sign (Q (z)) can admit by ana- 
lyzing the corresponding inequality systems. The idea is to look for sign patterns 
guaranteeing that every matrix with that sign pattern has a positive kernel vector. 
To this end we resort in subsection 13.21 to Qualitative Matrix Theory [7] and to L + - 
matrices in particular [23j . We first exemplify our approach by examining (|2.18l) for 
network (|2.1I) and turn to the general case in Section [ 
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3. Conditions for a singular reduced Jacobian G. 



3.1. System (|2.18p for network (f2TTj) . Note that for network (f2~Tj) the matrix 
E of (I2.8al) is also a basis for ker (S) (in general this need not be the case). Using this 
E we obtain for equation (|2.16[) (where gray indicates rows belonging to Q(z)): 



z 4 





















Z5 





























-Z6 













Z2 























Z6 




Z\ + 2z 3 - Z 4 











zi 4 


-2z 3 - 


" Zq 


z 2 - 


(- z 3 - z 5 








Z 2 - 


\-z 3 - 


Z6 














Zi - 


Ms - 


Z6 








Zl 







^6 













-22 




-Z6 





One has s = 5 and r = 10, hence the matrix Q(z) is defined by rows 6-10. However, 
it is easy to see that u£V (and hence positive v by (|2.8b[) ) exist only if z 6 = z 4 + z 5 . 
Hence Q(z) consists only of the rows 6, 7, 9 and 10 as row 8 evaluated at z 6 = z 4 + z 5 
is identically zero. One obtains 



Q(z) = 





zi + 2z 3 - z 4 













zi + 2z 3 — Zi — 







zz^t 


- z 3 - 


z 5 





Z2 + Z 3 — Z 4 — 












Zl 





z 4 + z 5 















-Z2 


-Zi - z 5 


one has positive 


iff the following 


I pairs 


of linear forms 


or both equal to 


zero: 










h(z) := Zl +2z a 


- z 4 


and 


h(z) 


= zi 4 


2z 3 - Zi- z 5 , 


£ 3 (z) := z 2 + z 3 


- z 5 


and 


h(z) 


= z 2 + 


Z 3 - Zi - z 5 , 




t 5 (z) : 


= Zl 


and 


e 6 (z) 


= z 4 4 


Z5, 




tr(z) := 


-z 2 


and 


h(z) 


= -z 4 


~ z 5 ■ 



Z5 
Z5 



(3.1) 



These conditions can be expressed as linear inequality systems, for example 

zi + 2z 3 - z 4 > 0, zi + 2z 3 — z 4 - z 5 < 0, 
Z2 + z 3 - z 5 < 0, z 2 + z 3 - z 4 - z 5 > 0, 
zi > 0, z 4 + z 5 < 0, 
— z 2 < 0, — z 4 — Z5 > 0. 



(3.2a) 



This system is feasible; pick any z £ M 5 satisfying (|3.2ap and let z§ — z 4 + Z5. Then 



Q(z) evaluated at that z has a positive kernel vector v (cf. Table l3~T1 and (|3.3[) . 
We apply Lemma E3] with M T = Z T = (0, 0, 1, 2, 1, 0) from (|2~3b| and need to find 
a vector Cj € ker (Z T ) with sign(cL)) = sign (z). For the choice of S) with uj 3 < 0, 
Cj 4 < and 0)5 > we consequently add 



z 3 < 0, z 4 < 0, z 5 > 
10 



(3.2b) 





1 


2 


3 


4 


5 


6 


z 


4 


1 


-5 


-8 


3 


-5 


UJ 


1 


1 


-3 


-1 


5 


-1 


X 


l 

4 


1 


3 
5 


i 

8 


5 
3 


i 

5 



Table 3.1 
Vectors z, ui and x 



to the inequality system. The overall system (|3.2ap & (|3.2b|) is feasible and one 
solution z is given in Table 13.11 Table 13.11 also contains a vector & G ker (Z T ^j with 
sign (z) — sign (cD) and the vector x — j (cf. Lemma [2~3l equation (|2.14|l ). Evaluating 
Q[z) at z from Table [3TT1 one has the matrix 



QO) 





-7 

4 





-1 

1 

-5 

-1 5 



that has the positive kernel vector 

v = (14, 4, 35, 140, 28) T . (3.3) 
Vector x from Table l3~T1 and the above v define a vector of rate constants: 



1400 160 12 672 

— , 112, — , — , — , 140, 63, 168, 140, 140 



Evaluation of Jac(k,x) at this k and x from Table [3~T1 confirms A = as a defective 
eigenvalue. 

All in all there are 81 different inequality systems where the pairs from (|3. 1[) are 
of different sign or both zero. There are also 13 inequality systems like (|3.2b[) that 
constrain z such that there is a u> <E ker (Z T ) with sign(w) = sign (z). Of these 
13*81=1053 inequality systems only the following four are feasible: 



z 3 < 0, z 4 < 0, z 5 > 
h{z) >0, £ 2 (z) <0 
£ 3 (z) < 0, t A {z) > 
£ 5 (z) > 0, £ e (z) < 
t 7 {z) < 0, £s(z) > 



z 3 < 0, z 4 > 0, z 5 < 
lx{z) < 0, t 2 {z) > 

(Pf ) £ 3 (z) > 0, U{z) < 

£ 5 (z) > 0, £ e (z) < 
£ 7 {z) < o, 40) > 



and 



z 3 > 0, z 4 > 0, z 5 < 
li{z)<Q, £ 2 {z)>0 
£ 3 (z) > 0, £i{z) < 

40) < o, 40) > o 
40) > o, 40) < o 



{Pi 



z 3 > 0, z 4 < 0, z 5 > 

40) > o, 40) < o 
40) < o, 40) > o 
40) < o, 40) > o 
40) > o, 40) < o 



(P 2 



n 



Because of the definitions of £5,^ and I7 in (|3.1j) . feasible z's do not have vanishing 
components. All in all we have established the following necessary and sufficient 
condition for a defective eigenvalue of J(v,x) of (12.11) . 

Fact 2. The Jacobian J(v,x) of \2. 1\) evaluated at (v,x) G V x -R§, (and hence 
Jac(k,x) evaluated at positive (k,x) via \2.6\) ) has X — as a defective eigenvalue, 
if and only if v and x satisfy: 

1. The vector x can be written as x — — with (i) z £ M 6 satisfies one of 
the inequality systems (P^), (P^) an ^ z 6 — z i + z 5 (implying Z\ ^ for 



i = 1, 



,6), (ii) lu G ker (Z T ^j and (Hi) sign (z) = sign (w). 



2. The above z and the vector v > are such that Q(z) v = 0. 

Note that, if z £ J£ 6 with zq = Z4 + Z5 satisfies one of the systems (P^) and (-P^), 
then the sign pattern sign(Q(z)) is one of the following: 
• If z satisfies (Pj 4 ) then 



sign(g(z)) = ± 



If z satisfies (P^) then 



sign(Q(.z)) = ± 



3.2. A sufficient condition. For the example of network (|2.1j) we obtained 
necessary and sufficient conditions in form of linear inequalities in z guaranteeing a 
positive kernel vector v of Q{z). The idea is to look for sign patterns sign(Q(z)) 
guaranteeing that every matrix with that sign pattern has a positive kernel vector (as 
it has been the case with the sign patterns of the previous section). By Qualitative 
Matrix Theory [7] (see in particular 23 ) one has the following Theorem 13.ll In our 
application, it can be stated in the following way: If a sign pattern sign(Q(z)) is an 
L + -matrix, then every matrix with the same sign pattern has a positive kernel vector. 

Theorem 3.1 (cf. [23], Theorem 2.4, p. 6). For a (m x n) sign pattern A, the 
following are equivalent: 

(a) A is an L + -matrix. 

(b ) Every matrix with the sign pattern A has a positive null vector and A has no zero 
row. 

(c) For each nonzero vector a G { — 1, 0, 1}™, some column of diag (a) A is nonzero 
and nonnegative. 

(d) For each nonzero vector a G { — 1, 0, l} m , some column of diag (a) A is nonzero 
and nonpositive. 

Note that Theorem 13.11 already contains - by the parts (c) or (d) - a primitive 
algorithm to determine whether or not a given sign pattern is an L + -matrix. So by 
Theorem 13 . 1 1 one can decide whether or not a particular sign pattern is an L + -matrix. 



With respect to the z-linear matrix Q(z) = Qij(z) G M sxp from (|2.17l) we propose 
the following: We first stack the columns of Q and consider the column vector 



(Q11, ■ ■ ■ , Qsi, 



■ j Qxpi ■ 



12 



Then we omit the components that are trivial linear forms to obtain a bijective map- 
ping of the form 

ip: Q(z)GM sxp h> Lz = {hz, J^zfeM 1 (3.4) 

with nontrivial n-dimensional row- vectors £i, . . . , £y, 7 < sp. So, the (7 x n)-matrix 

L = (4-) i= i,..., 7 ■ (3.5) 

just corresponds to the nontrivial linear forms in Q(z). Since we look for ai/£V with 
Q(z)v = 0, we are interested in the sign patterns that Lz can assume. So we define 
the set C + of all sign pattern matrices £ G { — l,0,l} sxp that are L + -matrices and 
that are realized by Q{z) for some z. Since the mapping (|3.4p associates a signature 
vector a = ^>(£) G {-1, 0, l} 7 to £ £ {-1, 0, l} sxp one arrives at 



C+ := {£ £ {-FO,!} 5 ^ 



£ is an L -matrix, 

(3.6) 

3z G Ft n with a % [Lz) i > if <7i ^ and (i z) 4 = if ct, = L 

Fact 3. Assume C + is nonempty and let £ G C + and a = V'(£)- Then there 
exists a vector z G J?" 

<7i (L z). > 0, i/ en ^ 0, (i z) . = 0, if<Ti=Q 

so that a — sign(Lz). Moreover, for each such z £ M n , there exists a positive 
v = v{z) with Q(z) v = by Theorem \3.1\ By the discussion of this implies 

that the pair (z, v) satisfies H(y) z = 0. 

Theorem 3.2. Consider a dynamical system defined by a mass action network 
as described in Section^ Recall the matrix Q(z) defined in H2.16}) . the matrix L 
defined in ^3.4^ & (375\) and the set C + defined in \3. 6\) . If there exist an element 
£ = sign(Q(z)) G C + and an element uj G ker (-Z T ) with 

sign(w) = sign(z), (3.7) 

then there exists a solution v G -K> t° Q( z ) v = and the Jacobian G(£, 77, v) of 
the reduced system has zero as an eigenvalue with algebraic multiplicity > 1. The 
corresponding steady state in original (k, x)- coordinates is given by 

x = ( x i)i=i,...,n> ( 3 - 8a ) 

Xi M' ifZi ^' (3.«») 
I Xi > ; arbitrary, if Zi — , 

fc = diag((/)(x- 1 )) Ev. (3.8c) 

The corresponding (£,f?)- coordinates in G{^,r],v) are then given by h2.10a\) . 
Proof. The statements follow directly from Lemma 12.31 and Fact [3J □ 



The condition (|3.7|) can be tested by examining the following linear inequality 
systems defined by orthants of R n . To establish these we identify each orthant by 
its sign pattern 5 G { — 1,0, 1}™: let x G M n , then the sign pattern of x is defined as 
5 := sign(a;) and the orthant containing x is given by Mg :— {x G IR n \ sign(x) = 5}. 
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To find z and ui satisfying (|3.7[) for a given signature a = ip(E) for an L + -matrix E 
then amounts to finding an orthant Mg such that 

<n (Lz) l > if <n ± 0, (Lz) l = if en = 0, (3.9a) 
Z T uj = 0, with <Jj cj, > 0, Si Zi > if S t f and ^ = 0,^=0 if <5; = 0. (3.9b) 

Corollary 3.3. If there exists a signature £ G C + and an orthant M$ , such 
that the linear inequality system 113. 9a)) , i3. 9b)) is feasible, then the reduced Jacobian 
G(£, v) has zero as an eigenvalue with algebraic multiplicity > 1. 

Remark 5. In the previous discussion we have only considered positive kernel 
vectors ofQ(z). However the setV can contain nonnegative vectors v. Thus, suppose 
v contains the facet of -R> given by \y G _ZR> |^ = 0, Vj > 0, i ^ j = 1, . . . , . 

Then one may fix v\ = ; replace Q(z) in the discussion above by the submatrix Q(z) 
obtained by deleting the i-th column (and eventually occurring zero rows) and obtain 
the remaining vi by asking for positive kernel vectors of Q{z) (i.e. by establishing the 
L + -property for Q(z) ). 

Remark 6. The condition \3. 9a)) tests whether the given L + -matrix S belongs 
to C + . By the definition of the matrix L this requires the labeling of the hyperplane 
arrangement given by Lz, which is computationally expensive (for an algorithm see 
Jl]j, J28}j). We have shown that all z satisfying liS. 9a)) lead to a positive null vector of 
Q(z). The condition \3. 9b)) then stands for the compatibility with the kernel of Z T : It 
tests whether there is a z in the solution set of 13. 9a)) that possesses a signature thats 
is compatible with ker (Z T ^j . Since one has to decide whether or not one of the systems 
\3.9a\ ). 13. 9b)) is feasible the overall procedure can be computationally demanding, even 
though the individual steps only involve simple matrix computations. 

4. Saddle node bifurcations for the Gl/S transition in budding yeast. 

The networks displayed in (|4.1I) and (|4.2[) below are competing hypotheses describ- 
ing the Gl/S transition in budding yeast. Both networks are biologically plausible 
and hard to distinguish experimentally |12j . 

[SiclP] - k3 -^- [0] [Sicl] 
[Clb] + [Sicl] - k " ' [Clb ■ Sicl] 




[Clb] + [SiclP] [Clb ■ SiclP] 

[Clb ■ Sicl] + [Clb] [Clb ■ Sicl ■ Clb] — H [Clb ■ SiclP] + [Clb] 

kn 

[SiclP] + [Cdcl4] =S^= [SiclP ■ Cdcl4] [Sicl] + [Cdcl4] 

ki4 

[Clb ■ SiclP] + [Cdcl4] [Clb ■ SiclP ■ Cdcl4] [Clb ■ Sicl] + [Cdcl4] 

ki7 
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[SiclP] [o] [Sicl] 



[Sicl ■ Clb] - k4 ~ [Clb] + [Sicl] - k " ~ [Clb • Sicl] 




kio 




(4.2) 



[Clb] + [SiclP] [Clb ■ SiclP] 

kn 

[SiclP] + [Cdcl4] [SiclP ■ Cdcl4] [Sicl] + [Cdcl4] 

ki4 

[Clb SiclP] + [Cdcl4] [Clb ■ SiclP ■ Cdcl4] — ^ [Clb ■ Sicl] + [Cdcl4] 



Switching is a desired property of models describing the Gl/S transition: depend- 
ing on its past a trajectory should move to different regions of state space, associated 
with the Gl and the S phase of cell cycle. Classically this has been realized by choosing 
rate constants and total concentrations, such that the ODE system shows bistability 
and hence hysteretic behaviour [51 [5Tj. For example in [12], multistationarity has 
been established for both models, indicating that both may be valid models. Here 
we demonstrate the applicability of our results by confirming switching for both net- 
works. We show that both models satisfy the conditions of Theorem [321 and compute 
states and rate constants where the Jacobian has a defective eigenvalue. We verify 
by numerical continuation that the system undergoes a saddle-node bifurcation, as 
generically expected, so that the codimension-1 stable manifold of the saddle-node 
and - after bifurcation - the one of the saddle represents a switching surface. 

For the network given in (|4.1[) one obtains using the stoichiometric matrix S given 
in Appendix |B| 



Z2-Z9 



-Z B + Z 

23 + 24 - 



Zi Z, + 2 3 - 2 9 



23 + 24 " 





23 + 24 
22 + 26 - 27 



3 —2g + 29 

29 23 + 24 — 29 





From the last three rows of Q(z) one has that positive v with Q(z) v = exist only if 



Z3 + Zi - Zg = 0, z 2 + z 6 - z 7 = 0, z 5 + z e - z 8 = 
and hence, for example, 



(4.3) 



Z-2 



-ZB + Z7, Z 3 



-Zi + Zg, Z 5 



-26 + Z8 



In this case colum 4, 5 and 6 will be the zero column, indicating that 1/4, 1/5, i>s > 
are unconstrained. Thus we need only consider the matrix 



Q s {z) 



-n$z 












-n w z 


-■Kg Z 


— 7T10 Z 


7T 6 Z 


7T 5 Z 








7T 2 Z 





— 7T4 Z 


7T 6 Z 


TT 5 Z 


7T 3 Z 





7T7 Z 


7T 3 Z 


7T7 Z 


7T7 Z 


7T 5 Z 





7T1 Z 


7T4 Z 
15 





7T4 Z 


IT 4 Z 


-7T 5 Z 



with 



7Tl 2 := — Z4 + 27 — 2 8 + 2g 
7T3 Z := Zl — 2Z4 + Zg 
7T5 ^ := Z 8 - Zg 
7T7 z := Zi ~ Z4 
7Tg 2 := Z4 



7T 2 Z := Z 5 + Z 7 - Z 8 - Zg 
7T4 Z != — Z5 + Zg 
7T 6 Z := Z 7 - Zg 
7T 8 Z := Z! 
7Tl 2 := Zg 



For example the system 



TV 1 Z < 0, 7T2 Z > 0, 7T3 Z < 0, 7T4 2 > 0, 7T5 Z < 0, 
7T6 Z > 0, 7T7 2 > 0, 7T 8 2 > 0, 7Tg Z > 0, 7Tlo 2 < 
2l > 0, 2 2 > 0, 2 3 < 0, 2 4 > 0, 2 5 < 0, 2 6 > 0, 2 7 > 0, 2 8 < 0, 2g < 

> 0, u>2 > 0, cj 3 < 0, lu 4 > 0, u; 5 < 0, loq > 0, uj-j > 0, oj 8 < 0, u;g < 



(4.4) 



and flU) hold. Then 



1 








1 


-1 


1 


1 


-1 











1 





-1 


1 


-1 





-1 





1 


-1 


1 


1 


-1 








-1 


1 





1 


1 


1 



is feasible. Let z E M 9 such that 



sign (Q s (z)) = 



is an L + -matrix (cf. (c,d) of Theorem 13. ip . One obtains, for example, the feasible 
points 

2 = (13, 2, -10, 8, -6, 2, 4, -4, -2) T , u = (1, 1, -1, 6, -1, 1, 1, -2, -1) T . 
Vectors z and Q yield the state vector 



x = — = 

z 



1113 11111 

13' 2' 10' 4' 6' 2' 4' 2' 2 



(4.5) 



For the matrix Q(z) evaluated at z one has 
Q(2) = 

The kernel of (5(5) contains the following positive vector 



-13 

















2 


-8 


2 


G 


-2 




















4 





-4 


G 


-2 





-5 














5 


-5 


5 


5 


-2 








-2 











4 





4 


4 


2 



18, 1, 2000, 1, 1, 1, 1. 



246 5114 1996 23146 Y 



35 105 



21 J 



Vectors x and v yield the rate constants 

k = 



1121 178204 4 328 72006 10228 1303760 

15 ' ' ' 3 ' 3' 35 ' ' 5 ' 35 ' 63 ' 



„ 65146 8004 , 7984 92668 „ 46292 \ 
2, — — , -— , 4, -— , —— , 2, 



(4.6) 



21 



21 
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A numerical continuation with this k and initial condition x verifies that at the dy- 
namical system undergoes a saddle-node bifurcation at (fc, x), cf. Fig. 14.1( a). 



For the network given in ()4.2p one obtains using the stoichiometric matrix S as 
given in Appendix [C] 



h Z3 — 24 
2 2 + 23 - 






1 + Za - 24 



2 2 - 2g 
2l - 23 + 29 





27 - z 9 

-21 -23 + Sg 




22 + 2 6 - 27 





-zg 
25 - 29 

-21 -23+29 



22 + 23 



-2 4 + 2 S - 2 9 

28 - 29 
-21 -23+2Q 



2S 22 + 2 3 - 28 



-21 -23+2g 

-24 + 2 8 
22 + 23 - 2 8 



From rows 3, 6 and 7 of Q(z) one has that positive v with Q(z) f = exist only if 

- zi - 23 + z 9 = 0, z 5 + z 6 - z s = 0, z 2 + z 6 - z 7 = (4.7) 

and hence, for example, 

z\ = -z 3 + z 9l z 2 = z 5 - z s + z 7 , z 6 = -25 + z 8 . 

In this case colum 2, 5 and 6 will be the zero column, indicating that z/2, ^5, ^6 > 
are unconstrained. Thus we need only consider the matrix 



Q s {z) 



7T 10 Z 








-7TH Z 


-7T 12 2 


7T 7 2 


-7T 12 2 


7T 4 2 


7T 9 2 














7T 3 2 


7T 7 2 


7T 8 2 


7T 9 2 


7T 9 2 





7T 5 2 





7T.5 2 











7T 6 2 


7T 6 2 








7T 2 2 











7T 2 2 


7Tl 2 


7Tl 2 



with 

7Ti 2 = 2 3 + 2 5 + 2 7 - 22 8 
7T3 2 = 2 5 + 2 7 - 2 8 - 2 9 
7T 5 2 = —24 + 2g 
7T7 2 = 27 — 29 
7T 9 2 = 2 8 - 2 9 
7Tn 2 = 24 

One obtains the feasible inequality system 



7T2 2 = 2 3 + 2 7 - 2 8 
7T4 2 = —24 + 2 8 — 2g 
7T 6 2 = -2 4 + 2 8 
7T 8 2 = 2 5 - 2g 

7T10 2 = 2 3 - 2g 

7T12 2 = 2g 



7Tl 2 > 0, 7T 2 2 < 0, 7T3 2 > 0, 7T4 2 = 0, 7T5 2 > 0, 7T6 2 < 0, 7T7 2 > 0, 
7T 8 2 > 0, 7Tg 2 < 0, 7Tl 2 < 0, 7Tn 2 < 0, 7Tl 2 2 < 0, 
2l > 0, 2 2 > 0, 2 3 < 0, 24 < 0, 2 5 > 0, 2 6 < 0, 27 > 0, 2 8 < 0, 2g < 0, 

u>i > 0, LU2 > 0, 103 < 0, LU4 < 0, > 0, lug < 0, W7 > 0, u; 8 < 0, cog < 0, 
where Q{z) is an L + -matrix. For example the feasible points 

z = (9, 9, -11, -4, 2, -8, 1, -6, -2) T , cD = (1, 1, -1, -1, 4, -1, 2, -1, -1) 
define the state vector 



(4.8) 



11 1 1 2 1 2 1 1 

9' 9' IT' 4' '8' ' 6' 2 
17 



(4.9) 



(a) Continuation for network 114,111 



(b) Continuation for network 1 14.211 



Fig. 4.1. Numerical continuation for the networks 114. H and (14.21 ) using rate constants k and 
initial condition x given in \j. 6\ l for network and in for network In 

both cases the upper and lower branches correspond to exponentially stable steady states. The total 
concentration ci is used as a bifurcation parameter. 



For the matrix Q(z) evaluated at z one obtains 

Q(z) = 



-9 



2 





0004 23 2 -4 

11 3 4 -4 -4 

0002 00 0-2-2 

4000 00 -4 4 4 



with the positive kernel vector 

v = (297, 11, 22, 11, 11, 11, 440, 1, 11, 451, 456, 6) T . 
Finally one obtains for the rate constants 

k = (1645, 2673, 9, 22, 92664, 45738, 112, 3584, 1850, 91476, 
Mi 1584, H ^, 1892, 66, 2772) T . 



(4.10) 



Again, numerical continuation show a saddle node bifurcation at (fc, x), cf. Fig |4.1f b). 
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Appendix A. Saddle-node bifurcations in mass action networks. 

This subsection is a recollection of some remarks concerning saddle-node bifurca- 
tions in mass action system that were originally made in [5] . We repeat them here to 
demonstrate the tight connection between zero eigenvalues of G(£, rj, A) obtained via 
(|2.12a[) fc (|2.12b[) and saddle node bifurcations. 

The following well-known theorem gives necessary and sufficient conditions for a 
saddle- node bifurcation of the system £ = g(£,r),k) as defined in (|2.10bj) . Let n be 
the component of the parameter vector v = (77, k) that will be used as bifurcation 
parameter. 

Theorem A.l (see e.g. [8], p. 497). Suppose (£0,^0) is a zero of g and suppose 
that the s x s matrix G(£o J l/ o) = S^Oi^o) has an algebraically simple eigenvalue 
with right eigenvector b and left eigenvector f3 T . Furthermore suppose that the 
following conditions are satisfied: 

g(£ ,v ) ^0, li T [Dlg{i ,u )(b,b)\ ^ 0. (A.l) 

Then there is a smooth curve of zeroes of g passing through (^cb^o)- Depending on 
the signs of the expressions in (|A.1[) , there are no or two zeroes near £0 for \i 7^ /Jq 
when the other components of vq remain fixed. 

In the remainder of this section we examine in terms of the Jacobian J(u, x) 
when the conditions (jA.lj) are satisfied for the reduced system (|2.10bj) . We recall the 
relations (|2.5aj) . (|2.6j) and (|2.10b|) . First observe that Fact 12.21 states conditions guar- 
anteeing that G(£, 7], A) has zero as an eigenvalue. For general mass action networks 
its algebraic multiplicity is expected to be 1. 

Now, if G(£o> Vo, v o) has an algebraically simple eigenvalue 0, we add two com- 
ments regarding (jA.ip (cf. [9]): 

(Nl) Recall the monomial function <fi(x), cf. Section [2] If fi is any rate constant 
fcj, then we have 

D n9(€o,fto) = U T S diag(0(x o )) 

and therefore 

fD^ gfa, no) = & (x ) /3 T U T Se^0 

for at least one i (as <j>i(xo) > and [U] = im(5)). 
(N2) From the above Lemma [2~T1 we deduce 

P T [Dlg(Zo,Vo,k o )(b,b)]^0^ 
a T [D 2 x f(x(Z , Vo ),ko)(A a ,Aa)] ^0 

with A :— J(A, x) and with left and right principal vectors a T and a of 
J(A,x). 

As a consequence of this discussion, in particular of the comment (Nl), we obtained 
in [9] the following remark concerning the originally given system (|2.2p , (|2.3al) : 

Remark 7. The system (|2.2|) . (|2.3aj) has a saddle-node bifurcation at (ko,xo) 
(within the plane Z T x = Z T xq =: c) if the following conditions are satisfied: 

(a) is a defective eigenvalue of J(Ao,a;o) with m a i g = m geo + 1 and the remaining 
eigenvalues have negative real parts. 

(b) a T \_D 2 X /(xo, ko)(Aa, Aaj\ ^ is satisfied for left and right principal vectors Wq 
and vq of J(Aq, xq). 
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Appendix B. The data for network (14.11) . 
B.l. Species and complexes of network (|4.1[) . 



Species 


Xi 


Complex 


Vi 


Sicl 


Xl 





2/1 


SiclP 


X 2 


Sicl 


2/2 


Clb 


X3 


SiclP 


2/3 


Clb ■ Sicl 


Xi 


Clb + Sicl 


2/4 


Clb ■ SiclP 


x 5 


Clb ■ Sicl 


2/5 


CdcU 


x 6 


Clb 


2/6 


SiclP ■ CdcU 


x 7 


Clb + SiclP 


2/7 


Clb ■ SiclP ■ CdcU 


x s 


Clb ■ SiclP 


2/8 


Clb ■ Sicl ■ Clb 


Xg 


Clb ■ Sicl + Clb 


2/9 






Clb ■ Sicl ■ Clb 


2/io 






Clb ■ SiclP + Clb 


2/ii 






SiclP + CdcU 


2/12 






SiclP ■ CdcU 


2/13 






Sicl + CdcU 


2/14 






Clb ■ SiclP + CdcU 


2/15 






Clb ■ SiclP ■ CdcU 


2/16 






Clb ■ Sicl + CdcU 


2/17 



B.2. Ordinary differential equat 



ions. 



X\ = k\ — ki xi — k 

X3 



: 4 Xi X 3 + k 5 X4 + ki5 X7 

-fc 3 x 2 - k 7 x 2 x 3 + k s x 5 - ki 3 x 2 x 6 + ku x 7 
-fc 4 xi x 3 + k 5 x A + k e Xi - k 7 x 2 x 3 + k 8 x 5 

+ kg X 5 - fcio £3 Xi + ku Xg + fci2 Xg 

ki xi X3 — k§ Xi — k$ Xi — kio X3 Xi + ku xg + fcig x% 
x 3 - k$ X 5 - kg X 5 + fci2 Xg - k m x 5 x 6 + k u x 8 
13 x 2 x 6 + ku x 7 + ki 5 x 7 - k 16 x 5 x 6 + k 17 x & + hs x 3 
kl3 x 2 x e - ku x 7 ~ ku x 7 
h 7 x 8 - km xs 

kll Xg — ki 2 Xg 



: k 7 X 2 



Xi — 

±5 = 

± 6 = 
X 7 = 

is = k 16 x 5 x 6 
Xg = fcio X3 Xi 



B.3. Conservation relations. 



Z{x 
Z 2 x 



- xe + x 7 + xs = ci 

-X3 + Xi + X 5 + X 8 + 2 Xg 



C'2 



B.4. The stoichiometric matrix: 



-llOOOOflflOO 



-1 -1 



1 -1 -1 -1 



-1 

1 



1 -1 
0-1 1 1-1 
1-1-1 




1 



1 -1 -1 

0000000001 -1 -1 0000 
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B. 5. The vector of reaction rates:. 

v(k,x) = (hi, k 2 x x , k 3 x 2 , k 4 xix 3 , k 5 X4, k 6 Xi, k 7 x 2 x 3l k$x 5l k 9 x 5 , k w x 3 Xi, 
kn%9, ki2Xg, k 13 x 2 x 6 , k 14 x 7 , k 15 x 7 , k 16 x 5 x 6 , k 17 x s , k 18 x s ) T 

Appendix C. The data for network ()4.2[) . 

C. l. Species and complexes of network (I4.2[) . 



Species 


Xi 


Complex 


Vi 


Sicl 


Xl 





yi 


SiclP 


x 2 


Sicl 


2/2 


Clb 


x 3 


SiclP 


2/3 


Clb ■ Sicl 


Xi 


Sicl ■ Clb 


2/4 


Clb ■ SiclP 


x 5 


Clb + Sicl 


2/5 


CdcU 


x e 


Clb ■ Sicl 


ye 


SiclP ■ CdcU 


x 7 


Clb 


V7 


Clb ■ SiclP ■ CdcU 


x 8 


Clb + SiclP 


ys 


Sicl ■ Clb 


Xg 


Clb ■ SiclP 


2/9 






SiclP + CdcU 


yw 






SiclP ■ CdcU 


yn 






Sicl + CdcU 


2/12 






Clb ■ SiclP + CdcU 


2/13 






Clb ■ SiclP ■ CdcU 


2/14 






Clb ■ Sicl + CdcU 


2/15 



C.2. Ordinary differential equations. 

xi = k\ — k 2 xi + fci xg — ks xi x 3 — ke xi x 3 + ^7X4 + ki§ x 7 
x 2 = -k 3 x 2 + kg xg - k w x 2 x 3 + fcn x 5 - ki 3 x 2 x 6 + ku x 7 
x 3 = &4 ig - &5 xi x 3 — kg xi x 3 + k 7 x 4 + ks £4 

+ k 9 x 9 - kw x 2 X 3 + fcn X 5 + kl2 x 5 
X4 = ks xi x 3 — k 7 Xi — ks Xi + kis x$ 

X 5 = fcio X 2 X 3 ~ fcn X 5 - fci2 x 5 ~ ki 6 X 5 X 6 + kyj x 8 

%6 = x 2 x 6 + ku x 7 + ki5 x 7 - ki 6 x 5 x 6 + kn x$ + ki S x 8 
±7 = fci 3 x 2 x e - ku X 7 - fci 5 x 7 
±8 = ^16 x 5 x e - kn x$ - k w x 8 

Xg = —ki Xg + k§ Xl X 3 — kg Xg 



C.3. Conservation relations. 

Zi x = x 6 + x 7 + x s = Cl 

X ~X 3 + Xi + X 5 + X$ + Xg = C 2 

C.4. The stoichiometric matrix:. 

"1-1 1-1-1 1 

00- 10000 
1-1-1 1 

000001- 1 

s= 







-1 1 
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-1 



C.5. The vector of reaction rates:. 

v(k,x) = (ki, k 2 xi, k 3 x 2 , k^xg, k 5 xix 3 , k 6 xix 3 , k 7 x^ k 8 x±, k 9 x 9 , k w x 2 x 3 , 
hix 5 , k 12 x 5 , k 13 x 2 x 6 , k 14 x 7l k 15 x r , k 16 x 5 x e , k 17 x 8 , k 18 x 8 ) 



Appendix D. Some Linear Algebra. 



D.l. The hypothesis im(S) = im(J(i/, x)) in (12.91) . We consider the factoriza- 
tion J = S [yV] T diag {x- 1 ) from l[2"3c |l &: (p3d] where V = diag (Ev) and diag (a;" 1 ) 
are diagonal matrices with positive entries. The equality im(S) = im(J) is thus equiv- 
alent to im(S) = im(S[yV] T ) since the invertible factor diag (a; -1 ) can be discarded. 

Given such (n x r)-matrices S and B = yV of rank (S) = s and rank (B) = 
rank(Y) =: f3 respectively, one always has im(S') D im(SB T ) and s > r&nk(SB T ). 
We discuss the equality 



S(M r ) = im (S) = im (SB T ) = S(im5 T ) 
which is obviously equivalent to s = ra,nk(SB T ) = rank(_B5' T ) and to 
dim (ker(SB T )) =n-s = dim (kcr(BS T )) . 



(D.la) 



(D.lb) 



Obviously, (|D.la[) necessitates f) > s. Moreover one has im(5) C im (SB T ) if and 



only if [im(5)] = ker{S T ) Dker(B T S) = [im (SB T 

BS T £, = Q S T £ 



0. 



and thus if and only if 

(D.lc) 



The elements £ of the (n— s)-dimensional subspace ker (S T ) satisfy (|D.lcjl . Therefore, 
rank (S) = ia,nk(SB T ) is equivalent to B \i m (s T ) being an injective map. This can be 
reformulated in terms of matrices S c and J5 C , whose columns form a basis of im (S T ) 
and im (B T ) respectively: rank (S) = rank (SB T ) is equivalent to B^ S c £ 1R TXS being 
of full column rank s. Finally, this fact leads with the help of matrices S c and y c , 
whose columns form a basis of im (S T ) and im (y T ) respectively, to the following 
characterization: 

im (S) = im (J(v, x)) & yj diag (Ev) S c £ M rxs has full column rank s . (D.ld) 



D.2. The reduced system Q(z)v = in (|2.17p . When discussing H(z) v = 
with H(z) given by (|2.7bj) . we have introduced matrices So and S c such that the 
columns of So and S c form a basis of ker(S') and im(S' T ) respectively. Furthermore we 
have chosen the Moore-Penrose inverse Sq and a particular matrix S# with S# So = 



It-, 



namely S^ a 



Theses choices have led to the equivalence of 

H(z) v — and Q(z)v — with the corresponding a will be given by P(z)v (cf. the 
set-up for PUS)) , ([2~T6f and IpTf]) ). 

When considering a different basis representation Sq = SqRo of ker(S') for a 
regular matrix Rq € jR( r - s ) x ( r ~ s ) and when working with the general form of S# 
given by 



ngen _ 

o # - 







'Ro 1 A " 




r9 # i 






. R T C _ 







(D.2a) 
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for regular matrices R c 6 M sxs and arbitrary matrices A £ M^ r S ) XS J one arrives at 
the equivalence of S a = diag (3^ T z) E u (cf. (|2.15|> ) to 



"i? - 1 A " 




'P(z) V 




a 


< 




Q{z) v_ 








(D.2b) 



for some matrix A of suitable dimensions. Hence, for general S^ n one obtains (in 
analogy to (|2 . 1 T[) ) the condition 

RlQ{z)v = Q. (D.2c) 

As R c is regular, one has that any pair (z, v)G JR" x V satisfying (|D.2cl) also satisfies 
(|2.17p and vice versa. Hence we conclude that (|2.17p is independent from the chosen 
bases for ker(S') and im(S T ). The corresponding a depends on the choice of Rq and 
A as ()D.2bl) shows. 

D.3. Ad Remark [4j By Appendix ID.ll one has im (5) = im ( J) if and only if 
one of the (s x s)-minors of diag (Ev)S c is nonzero for the chosen v from the kernel 
of Q(z). We like to add that these minors are polynomials in the components of Ev 
of order not greater than s. By Remark |4] following Theorem 12. 4[ the u's might be 
varied locally. Such a variation might be employed to establish (|D.ld[) . 
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